subroutine splitting_local_L_F
    use commondata
    implicit none
    integer :: i
    real*8 :: lamax
    do i = -4, im+4
        lamax = dmax1( abs(var(i,2)), abs(var(i,2)-var(i,4)), abs(var(i,2)+var(i,4)) )
        fp(i,1) = 0.5d0*( var(i,1)*var(i,2) + lamax*v1(i,1) )
        fp(i,2) = 0.5d0*( var(i,1)*var(i,2)**2 + var(i,3) + lamax*v1(i,2) )
        fp(i,3) = 0.5d0*( var(i,2)*( v1(i,3) + var(i,3) ) + lamax*v1(i,3) )
        fm(i,1) = 0.5d0*( var(i,1)*var(i,2) - lamax*v1(i,1) )
        fm(i,2) = 0.5d0*( var(i,1)*var(i,2)**2 + var(i,3) - lamax*v1(i,2) )
        fm(i,3) = 0.5d0*( var(i,2)*( v1(i,3) + var(i,3) ) - lamax*v1(i,3) )
    end do
end subroutine splitting_local_L_F


subroutine splitting_global_L_F
    use commondata
    implicit none
    integer :: i
    real*8 :: lamax
    lamax = 0.d0
    do i = -4, im+4
        lamax = dmax1( lamax, abs(var(i,2)), abs(var(i,2)-var(i,4)), abs(var(i,2)+var(i,4)) )
    end do
    do i = -4, im+4
        fp(i,1) = 0.5d0*( var(i,1)*var(i,2) + lamax*v1(i,1) )
        fp(i,2) = 0.5d0*( var(i,1)*var(i,2)**2 + var(i,3) + lamax*v1(i,2) )
        fp(i,3) = 0.5d0*( var(i,2)*( v1(i,3) + var(i,3) ) + lamax*v1(i,3) )
        fm(i,1) = 0.5d0*( var(i,1)*var(i,2) - lamax*v1(i,1) )
        fm(i,2) = 0.5d0*( var(i,1)*var(i,2)**2 + var(i,3) - lamax*v1(i,2) )
        fm(i,3) = 0.5d0*( var(i,2)*( v1(i,3) + var(i,3) ) - lamax*v1(i,3) )
    end do
end subroutine splitting_global_L_F


subroutine splitting_Steger_Warming
    use commondata
    implicit none
    integer :: i
    real*8 :: la(3), la_p(3), la_m(3), wp, wm, temp
    do i = -4, im+4
        la(1) = var(i,2)
        la(2) = var(i,2) - var(i,4)
        la(3) = var(i,2) + var(i,4)
        la_p(:) = 0.5d0*( la(:) + abs(la(:)) )
        la_m(:) = 0.5d0*( la(:) - abs(la(:)) )
        temp = 0.5d0*var(i,1)/gamma
        
        fp(i,1) = temp*( 2.0d0*(gamma-1.0d0)*la_p(1) + la_p(2) + la_p(3) )
        fp(i,2) = temp*( 2.0d0*(gamma-1.0d0)*la_p(1)*la(1) + la_p(2)*la(2) + la_p(3)*la(3) )
        wp = 0.5d0*(3.0d0-gamma)*( la_p(2) + la_p(3) )*var(i,4)**2/(gamma-1.0d0)
        fp(i,3) = temp*( (gamma-1.0d0)*la_p(1)*la(1)**2 + 0.5d0*la_p(2)*la(2)**2 + 0.5d0*la_p(3)*la(3)**2 + wp )
        
        fm(i,1) = temp*( 2.0d0*(gamma-1.0d0)*la_m(1) + la_m(2) + la_m(3) )
        fm(i,2) = temp*( 2.0d0*(gamma-1.0d0)*la_m(1)*la(1) + la_m(2)*la(2) + la_m(3)*la(3) )
        wm = 0.5d0*(3.0d0-gamma)*(la_m(2)+la_m(3))*var(i,4)**2/(gamma-1.0d0)
        fm(i,3) = temp*( (gamma-1.0d0)*la_m(1)*la(1)**2 + 0.5d0*la_m(2)*la(2)**2 + 0.5d0*la_m(3)*la(3)**2 + wm )
    end do
end subroutine splitting_Steger_Warming


subroutine splitting_VanLeer
    use commondata
    implicit none
    integer :: i
    real*8 :: ma, mp, mm
    do i = -4, im+4
        ma = var(i,2)/var(i,4)
        if( ma >= 1.0d0 ) then
            fp(i,1) = var(i,1)*var(i,2)
            fp(i,2) = var(i,1)*var(i,2)**2 + var(i,3)
            fp(i,3) = ( v1(i,3) + var(i,3) )*var(i,2)
            fm(i,:) = 0.d0
        else if( ma <= -1.0d0 ) then
            fm(i,1) = var(i,1)*var(i,2)
            fm(i,2) = var(i,1)*var(i,2)**2 + var(i,3)
            fm(i,3) = ( v1(i,3) + var(i,3) )*var(i,2)
            fp(i,:) = 0.d0  
        else
            mp = 0.25d0*(ma+1.d0)**2*var(i,1)*var(i,4)
            fp(i,1) = mp
            fp(i,2) = mp*( (gamma-1.d0)*var(i,2) + 2.d0*var(i,4) )/gamma
            fp(i,3) = 0.5d0*mp*( (gamma-1.d0)*var(i,2) + 2.d0*var(i,4) )**2/(gamma**2-1.d0)
            mm = -0.25d0*(ma-1.d0)**2*var(i,1)*var(i,4)
            fm(i,1) = mm
            fm(i,2) = mm*( (gamma-1.d0)*var(i,2) - 2.d0*var(i,4) )/gamma
            fm(i,3) = 0.5d0*mm*( (gamma-1.d0)*var(i,2) - 2.d0*var(i,4) )**2/(gamma**2-1.d0)
        end if
    end do
end subroutine splitting_VanLeer